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We discuss the fluctuation properties of equilibrium chaotic systems with constraints such as iso- 
kinetic and Nose-Hoover thermostats. Although the dynamics of these systems does not typically 
preserve phase-space volumes, the average phase-space contraction rate vanishes, so that the sta- 
tionary states are smooth. Nevertheless finite-time averages of the phase-space contraction rate have 
non-trivial fluctuations which we show satisfy a simple version of the Gallavotti-Cohen fluctuation 
theorem, complementary to the usual fluctuation theorem for non-equilibrium stationary states, and 
appropriate to constrained equilibrium states. Moreover we show these fluctuations are distributed 
according to a Gaussian curve for long-enough times. Three different systems are considered here, 
namely (i) a fluid composed of particles interacting with Lennard- Jones potentials; (ii) a harmonic 
oscillator with Nose-Hoover thermostatting; (iii) a simple hyperbolic two-dimensional map. 
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Recently some papers have appeared which discuss 
the applicability of the Gallavotti-Cohen Fluctuation 
Theorem (FT) [l| to near equilibrium and equilibrium 
chaotic systems subject to isokinetic and Nose-Hoover 
thermostats 

H H II In particular, in 0, Evans et 
al. consider constrained isokinetic systems of interact- 
ing particles which are in equilibrium or near equilibrium 
states and claim the fluctutations of the phase-space con- 
traction rate of the systems at or near equilibrium do not 
satisfy the FT. 

One of the issues under discussion in these papers is 
whether or not the FT applies when the external forces 
vanish, which we refer to as constrained equilibrium. 
These systems differ from genuine equilibrium ones be- 
cause of the action of a constraint, for instance acting so 
as to keep the kinetic energy constant. For systems con- 
sisting of particles with instantaneous, elastic, hard ball 
collisions the total energy is all kinetic and the constraint 
has no effect in the absence of external forces, thus the 
resulting ensemble is the usual microcanonical one. For 
other systems the kinetic and total energies are not the 
same. The constraint fixes only the kinetic energy and 
the total energy fluctuates. As a result, the constrained 
equilibrium ensemble is not microcanonical, but has a 
different structure, which we explore here. 

For the sake of illustration, let us consider a fluid com- 
posed of particles with pair-wise, central interactions, 
placed in a periodic box, and not subjected to external 
forces. In the absence of other constraints, the micro- 
canonical measure is the invariant one. Now we mod- 
ify the dynamics by introducing a (time-reversible) fric- 
tional force on the particles with the requirement that 
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the total kinetic energy remain constant in time, see 
Eqs. ©, i|10|l an d |Ht below for a specific example. 
Because of this constraint, the total energy is not con- 
served and the Liouville theorem is violated. In such a 
case the phase-space flow can have non-vanishing diver- 
gence (point-wise), even though it vanishes in average. 
In other words, as long as there is no external driving, 
the system relaxes to an equilibrium state, different from 
the one specified by the volume measure, but neverthe- 
less with vanishing average phase-space contraction rate. 
For Anosov systems or those satisfying the Chaotic Hy- 
pothesis, this implies the sum of the Lyapunov exponents 
is identically zero. This is to say the stationary mea- 
sure is smooth both with respect to stable and unstable 
manifolds Moreover, the existence of a time-reversal 
symmetry implies that the equilibrium measure must be 
symmetric under time-reversal operation, i. e. a trajec- 
tory and its time reverse are equally as probable. This 
is indeed what one expects from an equilibrium state, 
whether it verifies Liouville's theorem or not. 

The aim of this paper is to provide through some el- 
ementary considerations a characterization of the fluc- 
tuations of the phase-space contraction rate of such con- 
strained equilibrium systems, which we will subsequently 
support through numerical studies of specific examples. 
The outcome indicates that the FT does apply, in agree- 
ment with conclusions reached by Gallavotti et al. 

For non-equilibrium stationary states, the FT is a 
statement about the asymmetric part of the probability 
distribution of the phase-space contraction rate. Let a = 
(a) denote the non-zero expectation value of the phase- 
space contraction rate of the non-equilibrium stationary 
state of a given system and let Prob(er t = pa) denote the 
probability of observing, during a time interval of length 
t, a given fluctuation of the partial average of a equal 
to pa, where the amplitude p is a dimensionless number. 
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The FT states that when a > 



1 Prob(cr t = per) _ 
lira — in — — — — = pa. 

t^oo t Prob(cr t = — pa) 



(1) 



As pointed out in |j] , the derivation of this result assumes 
W > and must be reformulated if W = 0. 

For an equilibrium system, finite-time distributions of 
the phase-space contraction rate are symmetric about 
zero. Now let a denote the point-wise divergence of a 
flow with the properties described above (or equivalently 
the logarithm of the Jacobian of a time-discrete mapping 
with similar properties), and let Prob(er = a) denote the 
(equilibrium) probability of observing a given value a of 
a. For an equilibrium system, the most general FT is the 
following : 



Prob(er 



Prob(cr = —a), 



(2) 



which holds since the equilibrium distribution is symmet- 
ric under time-reversal (this operation changes the sign 
of a). Equation (J2J is similar to detailed balance for 
stochastic systems. 

For a system whose phase-space volumes are preserved 
by the evolution, the probabilities in Eq. © are non- 
vanishing only at a — 0, so that no fluctuations occur. 
But this is not necessarily so for constrained equilibrium 
systems. If so, one can consider the ratio of the two 
probabilities, equal to unity, and average a over some 
time interval without affecting the result. Let a t denote 
the average of a over a time interval of length t and take 
the logarithm of the ratio of Eq. J5J) applied to er t to 
obtain the expression 



a) 



1 Prob(cr t 
lim — In — 

t^oo t Prob(<7t = —a) 



0. 



(3) 



which is a corollary to Eq. and, as remarked in £|, 
is the natural generalization of the Gallavotti-Cohen FT, 
here applied to an equilibrium system. The statement of 
Eq. is weaker than Eq. @ and not very useful for its 
own sake since the equilibrium stationary state carries 
no asymmetric part. Nevertheless it is correct and not in 
contradiction with the FT as stated in . 

This observation does not support the comments by 
Evans et al. Q that for Anosov equilibrium dynamics the 
range of the admissible fluctuations for the FT shrinks 
to zero. While it is true that the asymptotic range of 
admissible fluctuations of at may shrink to zero as t — > 
oo, one should still expect to measure fluctuations about 
the value a = as long as t < oo. Moreover the rate 
at which Prob(<r t = a) goes to zero is the same as that 
of Prob(<7t = —a), so that the ratio of the two remains 
constant as t — > oo. Of course, the right hand side of 
Eq. being zero, the result is trivial, irrespective of 
the asymptotic form of the probability distribution. In 
other words, unlike Eq. Q for non-equilibrium stationary 
states, neither Eq. ©, nor Eq. © hold any information 
on the asymptotic form of the probability distribution 



of the phase-space contraction rate; the only information 
provided by Eq. @ is that the probability distributions 
are symmetric about their average value, zero. Moreover, 
Eq. (2J holds for equilibrium systems irrespective of the 
details of the dynamics, whether it verifies the Chaotic 
Hypothesis or not. The details of the fluctuations may 
indeed vary from one system to another, but even though, 
whether the fluctuations are trivial or not, Eq. @ always 
holds. 

Here we will show, through the example of an equi- 
librium hyperbolic map [Toj with non-trivial Jacobian, 
that the distribution of the phase-space contraction rate 
is symmetric about zero, with an asymptotic form given 
by Prob(<7t = x) « pt(x)dx, with 



Pt(x) 



: exp 



{xtf 



2x 2 



(4) 



where \ is a constant which depends on the specific sys- 
tem under consideration. We note that this indicates 
that the distribution of the total phase space contraction 
accumulated in time has a variance that is asymptotically 
independent of time. 

We will show through two additional examples that 
for both Nose-Hoover thermostatted equilibrium systems 
(whose kinetic energy distribution is given by the canon- 
ical ensemble distribution) and constrained, equilibrium 
interacting particle systems, such as a Lennard- Jones 
fluid with isokinetic thermostat, the asymptotic form of 
the phase-space contraction rate distribution function is 
(i) non-trivial, (ii) consistent with the relevant Equilib- 
rium Fluctuation Relations, Eqs. (ll'IMII and (iii) verifies 
a central limit theorem similar to Eq. Q, albeit with a 
time-dependence in \fi instead of t. 

The first example we consider is a mapping of the unit 
square into itself with hyperbolic properties similar to the 
baker map. It is a variant of a non-linear baker map pre- 
viously introduced by us (3- Namely let M a be defined 
by 

M a (x,y)= (5) 

(2(p a (x),<p- a (y/2)), 0<x<l/2, 
(2<p a {x) - 1, <p- a ((y + l)/2)), 1/2 < x < 1, 



where 



<Pa(x) 



(6) 



(l/7r)arctan [tan(7rx)2~ a ] , < x < 1/2, 
1 + (l/7r)arctan [tan(7ra;)2~ Q l, 1/2 < x < 1. 



The time-reversibility of M a under the time-reversal 
operator T(x, y) = (1 — y, 1 — x) is easy to check. More- 
over, as long as < a < 1, M a is uniformly expanding 
along the x-coordinate, while it is contracting along the 
y-coordinate. It is smooth, i. e. it is continuous and 
has continuous derivative everywhere, but for the cut at 
x = 1/2. Therefore it shares the hyperbolic properties 
of the usual baker map. The Jacobian, which we de- 
note \DM a (x,y)\, is typically not equal to unity, but the 
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natural invariant measure of this map turns out to be ab- 
solutely continuous, as far as one can tell from numerical 
simulations. This is illustrated in Fig^where the param- 
eter was set to a — 1/2. The symmetry of the invariant 
measure under the time-reversal operator is clearly seen 
from the figure. 

The smoothness of the invariant measure is consistent 
with the computation of the Lyapunov exponents, which, 
for this value of the parameter, are numerically found to 
be A+ = 0.6472 + o(l(T 4 ), A_ = -0.6472 + o(l(T 4 ), i. e. 
within numerical errors A + + A_ = 0. Figure [21 shows 
a numerical computation of the probability distribution 
function of the logarithm of the Jacobian of M a , which 
is symmetric about 0, as expected from Eq. J5J. The 
probability distributions of the time-averaged Jacobian 
are shown in Fig. 1 3141 

The second example we consider is that of a harmonic 
oscillator with a so-called Nose-Hoover chain of ther- 
mostats here limited to two thermostat variables. This 
generalization of the Nose-Hoover thermostatting scheme 
is due to Q The system is time-reversible and yields 
canonical distributions of the phase variables, as shown 
in 0. 

Let q denote the position variable, p the momentum, 
£ and £ the thermostatting variables. Taking all the os- 
cillator parameters to be unity, the equations of motion 
are (see p. 193 of 0) 



g = p, 

p = -q- (p, 
C = p 2 -1-£C, 



(7) 



In the absence of a thermostatting mechanism, oscilla- 
tions would be periodic of period 2tt, but this is no longer 
so if the thermostatting mechanism is turned on. 

The phase-space contraction rate is a = — C — £ whose 
distribution, computed numerically with an adaptative 
step size 5th order Runge-Kutta scheme 0, is shown in 
Fig. This curve is Gaussian with standard deviation 
a/2, consistent with the analytical result given in 0, P- 
193. Distributions of time-averaged contraction rates at 
are shown in Figs. 16171 Deviations from an exact Gaus- 
sian can be seen for intermediate times, which are due to 
surviving correlations. However, these deviations disap- 
pear for long enough times, so that one retrieves a Gaus- 
sian distribution whose asymptotic form is now given by 
Prob((7t = x) sa pt(x)dx, with 



Pt{x) 



2ttx ; 



■ exp 



{x)H 



2x 2 



(8) 



The third example we consider is that of a Lennard- 
Jones fluid with the kinetic energy constrained to a con- 
stant value, 



1i = 

i>i = 



F(qij) - apt 



(9) 
(10) 



where F(q^) denotes the force between particles i and j 
associated to the usual Lennard-Jones potential, and 
is the distance separating the particles, 



F(r) = -V0(r), 

W = 4e [G) 12 -(D f 

Qij = |q» — qjl- 



(ii) 

(12) 
(13) 



The damping term a in Eq. IjlOfl is a function of the 
phase-space coordinates, derived from the condition that 
the kinetic energy be constant and such that the equa- 
tions I|9I10(I are time-reversal symmetric, 



Eij,« Pi-F(«y) 



(14) 



Notice the time dependence is different from Eq. @. 



One readily sees a is precisely minus the divergence of 
the flow (up to a multiplicative constant), and therefore 
the quantity whose probability distribution is expected to 
be symmetric about zero according to Eq. J2J. This has 
been checked numerically for a two-dimensional fluid of 
8 particles in a box of approximately 7x8 units of length 
squared, using the modified Verlet leapfrog algorithm, as 
described in [9j. 

The resulting probability distribution of a is shown in 
Fig.[S]and the probability distributions of its time average 
in Figs. lOTTOl 

In summary, we have shown for three examples of con- 
strained equilibrium chaotic systems, that the probability 
distributions of the phase-space contraction are symmet- 
ric about zero and verify a central limit theorem. The 
time-dependence varies according to the nature of the 
dynamics. It would be interesting to further investigate 
this property. 
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FIG. 1: (Color online) Numerical computation of the invari- 
ant density of M a {x,y) Eq. (JSJ for a — 1/2. The computation 
is the result of 10 4 trajectories iterated over a time of 10 3 
steps. The same time was discarded prior to outputting tra- 
jectories in order to eliminate transient effects. 




FIG. 2: Numerical computation of the probability density of 
the logarithm of the Jacobian of M a {x, y) Eq. © for a = 1/2. 
The computation is the result of sampling 10 4 trajectories 
iterated over a time of 10 4 steps. 10 3 bins were used. 
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FIG. 3: Numerical computation of the probability density 
of the logarithm of the time-averaged Jacobian of M a {x,y) 
Eq. (jJJ for a = 1/2. The curves shown here are averaged 
respectively over N = 5, 10, 15, 20, 25, 50, 75, 100, 150, 200, 
250, and 500 steps (the curves get narrower as the number of 
steps increases). The computation is the result of sampling 
10 4 trajectories iterated over a time of 10 4 steps. 10 3 bins 
were used. 



FIG. 4: Same as Fig. [3] with the curves rescaled according to 
Eq. J3J . The thick dashed curved is a Gaussian with standard 
deviation \ 2 — 0.35. 
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FIG. 5: Numerical computation of the probability density of FIG. 6: Numerical computation of the probability density of 



— £— £, the divergence of the Nose-Hoover oscillator, Eq. $Jj, 
for the trajectory with initial conditions are go = 0, po = 1, 
("o = 0.5, £o = 0. The integration is performed over a time 
2iy x 10 s , sampled every 2tv, with 500 bins spanning [—10, 10]. 
The thick line on top is a Gaussian with standard deviation 
y/2. 



the time-averaged contraction rate at, for times T = 1, 5, 10, 
15, 20, 25, 50, 75, 100, 150, 200, 250, 500 and 750 (times 2tt). 
The parameters of the integration are the same as Fig. 
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FIG. 7: Same as Fig. |S| with the curves rescaled according to 
Eq. ®. The thick curve is a Gaussian with standard devia- 
tion x 2 — 3-2- The agreement is excellent for times 200 and 
greater. 
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FIG. 8: Numerical computation of the probability distribu- 
tion of a, the damping term in Eq. llOH . The parameters of 
the simulation correspond to the usual choice of dimensionless 
parameters. 8 particles were used on a periodic 2-dimensional 
box of size approximately 7x8. The temperature was set to 
T = 1 and the particle masses to m = 1. A number of 
10 4 trajectories were taken and the value of a recorded for 
1000 units of time at intervals of 1 time units. 300 bins were 
used, spanning values a £ [—3,3]. The average value of a 
was computed to be « 2.7 x 10 -4 with standard deviation 
ss 1.8 x 10~ 3 . About 8.7% of the computed values fell out of 
the range shown in the figure. 
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FIG. 9: Same as Fig. [S]for the time-average at = 1/T J a, FIG. 10: Same as Fig. [5] with the curves rescaled according 
for times T = 1, 2, 3, 4, 5, 10, 15, 20, 25, 50, 75, 100, 150, to Eq. ©. The thick curve is a Gaussian with standard devi- 
200, 250, 500 and 750. ation \ 2 — 5.9. The agreement is excellent for times 200 and 

greater. 



